Critical behaviour of a tumor growth model — 
Directed Percolation with a mean-field flavour 
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We examine the critical behaviour of a lattice model of tumor growth where supplied nutrients 
are correlated with the distribution of tumor cells. Our results support the previous report (Ferreira 
et al., Phys. Rev. E 85, 010901 (2012)), which suggested that the critical behaviour of the model 
differs from the expected Directed Percolation (DP) universality class. Surprisingly, only some 
of the critical exponents (/?, a, v±, and z) take non-DP values while some others (/?', uu, and 
spreading-dynamics exponents O, S, z') remain very close to their DP counterparts. The obtained 
exponents satisfy the scaling relations f3 — avn, /?' = Svu, and the generalized hyperscaling relation 
Q + a + 8 — d/z, where the dynamical exponent z is, however, used instead of the spreading 
exponent z . Both in d = 1 and d — 2 versions of our model, the exponent /3 most likely takes the 
mean-field value /3 = 1, and we speculate that it might be due to the roulette- wheel selection, which 
is used to choose the site to supply a nutrient. 
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I. INTRODUCTION 

Recently, critical behaviour of models with absorbing 
states is attracting a lot of interest [3, [1] . A lot of such 
models appeared mainly due to a large number of phe- 
nomena where they are potentially applicable, such as, 
for example, in spreading of fires or epidemics, heteroge- 
neous catalysis or granular flows. Yet another motivation 
comes from the intriguing idea that critical behaviour of 
these noncquilibirum systems can be, at least to some 
extent, categorized and understood in an analogous way 
to equilibrium systems, where the notions of scale invari- 
ance, universality, and renormalization group turned out 
to be so fruitful and productive [3|. 

A prominent and strong indication that the universal- 
ity is indeed common also among nonequilibrium systems 
comes with Directed Percolation (DP). Defined in a sim- 
ple mathematical way, this spreading process turns out 
to exhibit the critical behaviour, which is very robust and 
appears in a number of other models as well as in some 
experimental realizations Q. A common feature of sys- 
tems exibiting the DP criticality is that they have a sin- 
gle (or effectively single) absorbing state. The ubiquity 
of the DP critical behaviour prompted three decades ago 
Janssen [j| and Grassberger Q to formulate a conjecture 
that all continuous transitions in systems with a single 
absorbing state and lacking extra symmetry or conserva- 
tion law will belong to the DP universality class. Up to 
now, there are plenty of numerical examples supporting 
this conjecture [l|, H|. Models with multiple absorbing 
states might also belong to the DP universality class M] 
but in some cases they exhibit a different behaviour Q. 

Absorbing states appear also in some ecological or tu- 
mor growth models, where they correspond to the ex- 
tinction of one or more species or types of cells. In some 
cases the critical points in these models were shown to be- 



long to the Directed Percolation universality class [9l-fl2j . 
However, an implementation of biologically relevant fac- 
tors, such as, for example, inter- or intra-species competi- 
tion, food supply, aging, birth or death, results in rather 
complex dynamics. Hence, a violation of the Grassberger 
and Janssen conjecture cannot be excluded in such mod- 
els. 

In the present paper, we examine the critical behaviour 
of a recently introduced tumor growth model 0. Our 
previous report [l3[ suggested that despite a single ab- 
sorbing state, the critical behaviour of this model dif- 
fers from the expected DP universality class. Having 
performed much more extensive and detailed analysis, 
in the present paper we confirm the non-DP criticality 
of the model. However, the apparent violation of the 
Grassberger and Janssen conjecture occurs in a very in- 
tricate way: only some exponents take non-DP values 
while some others, including all spreading-dynamics ex- 
ponents, coincide with their DP counterparts. Neverthe- 
less, some scaling and hyper-scaling relations seem to be 
satisfied in our model. 



II. MODEL AND NUMERICAL METHODS 

A. Definition of the model 

In our model, each site of a c?-dimensional lattice of 
linear size N either is occupied by tumor, nutrient, both 
tumor and nutrient, or is empty. At a rate p, nutri- 
ents are supplied to a chosen site of the lattice, provided 
that the site is not already occupied by a nutrient. The 
roulette- wheel selection [14J is used to choose the site for 
such a supply and the corresponding weight w depends 
on whether the site is occupied by tumor (w = 1 + A) or 
not (w — 1). The parameter A > takes into account 
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the angiogenic effects of increased nutrients supply due 
to formation of new blood vessels in the vicinity of tu- 
mor cells. At a rate 1 — p, a tumor cell is updated on a 
randomly chosen site. The tumor cell survives provided 
that there is a nutrient on this site, otherwise it dies. The 
surviving tumor cell consumes the nutrient and attempts 
to breed provided that there is a site without a tumor 
cell among its nearest neighbours. Some preliminary re- 
sults on this model were already reported Closely 
related models but without any preference in nutrient 
supply (A = 0) were also studied [13 ] . 



B. Monte Carlo simulatons 

To examine our model, we used Monte Carlo simula- 
tions. To facilitate more extensive simulations, we imple- 
mented an algorithm where a site is selected only from 
among active sites, i.e., sites that are either without nu- 
trients or with tumor cells. Maintaining the lists of such 
sites requires some bookkeeping but it speeds up the sim- 
ulations substantially, especially in the vicinity of the 
critical point (where most sites are inactive) . We checked 
that the simpler algoritm that does not keep the list of 
active sites and the present algorithm yield statistically 
the same results with respect both to static as well as to 
dynamic characteristics. 

For various p and A, we measured the steady-state 
densities of tumor cells Xt and of nutrients x n . Simu- 
lations started from a random initial configuration and 
the model relaxed until a steady state was reached. To 
ensure that the obtained results are N- independent, we 
examined lattices of various sizes. Assuming that xt is 
an order parameter, we expected that upon approaching 
the critical point (p = p c ), it decays as xt ~ (p — p c )^ , 
where j3 is a characteristic exponent. We also measured 
the time dependence of the tumor cell density xt (i) (with 
the unit of time corresponding to N d update attempts). 
The density Xt{t) for each t is an average over indepen- 
dent runs. At the critical point, the density x t (t) is ex- 
pected to have a power- law decay x* (t) ~ t~ a , where a 
is another characteristic exponent [l|, [l5| . 

Some other exponents can be extracted assuming the 
following scaling form of the tumor cell density [l[ 

x t {t)~t- a fUp-p c )t 1 '^,t d ''/N) (1) 

In particular, from the finite-size scaling analysis of the 
density x t (t) at the critical point (p = p c ), we can esti- 
mate the dynamical exponent z from the collapse of the 
suitably rescaled x t (t) for various N. With a similar pro- 
cedure, the exponent v\\ which describes the divergence 
of the temporal correlation length (£n ~ (p — p c ) _iy n), 
can be estimated from the off-critical steady-state data. 

The divergence of the spatial correlation length £j_ is 
characterized with the exponent vy ~ (p — Pc)~ l ' ± )- 
For finite-size systems at criticality this exponent governs 
also the vanishing of the steady-state order parameter 
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x t ~ N~^ v± (2) 

Further information about the critical behaviour of our 
model can be obtained using the so-called spreading dy- 
namics. In this technique one starts the simulation from 
a nearly-absorbing configuration and monitors the char- 
acteristics of the runs. In the case of our model, we 
started the simulations from the configuration where one 
site contains both a nutrient and a tumor cell and all 
other sites are filled only with nutrients. We measured 
the average number (not density) of tumor cells {Xt), the 
survival probability Pit) that a given run survives until 
time t, and the mean square of spreading of tumor cells 
from the origin R 2 {t). At criticality one expects that 

X t {t)~t e , P{t)~t 5 , R 2 {t)~t 2/z ' (3) 

where 0, S and z' are characteristic exponents of the 
spreading process. Moreover, from the decay of the long- 
time limit of the survival probability 

P(t oo) ~ (p - p c f (4) 

we can extract yet another characteristic exponent 

There are a number of examples showing that for mod- 
els with a single absorbing state, the dynamical exponent 
z and the spreading exponent z' are equal. A violation 
of z = z' equality was reported in some nonequilibrium 
Ising models but it was attributed to the missing ab- 
sorbing phase [l7j . Moreover, for models with multiple 
absorbing states, the exponent z' might depend on the 
absorbing state that is chosen in spreading-dynamics sim- 
ulations [isl | . 

The above critical exponents are not entirely indepen- 
dent. For models like Directed Percolation or Contact 
Process a — 5, z' = z, ft = j3' and some plausible scaling 
arguments can be used to derive the following hyperscal- 
ing relation [l9[ 

9 + 28 = - (5) 

z 

For models with multiple absorbing states the hyperscal- 
in;? relation is expected to hold in a more general version 

6 + a + 6 = 4 (6) 
z' 

and in general f3 ^ /3' . 

Using scaling assumptions one can also show that [l| 

/3 = oi/||, P' = 5v\\. (7) 

Moreover, the dynamical exponent z, which is defined 
from the relation £|| ~ must of course satisfy the 
equality 

z = vu/v x . (8) 
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C. Mean-field approximation 

To supplement Monte Carlo simulations, we investi- 
gated the behaviour of our model using the Mean-Field 
Approximation. The general strategy to develop such an 
approach is to write a master equation and then solve 
it using some decoupling procedure [H, US [HI • in the 
first step one sums the states of all but one site. As a 
result, one obtains equations describing the time evolu- 
tion of probabilities for single-site configurations. Such 
equations will contain probabilities of more complex con- 
figurations (i.e., two-site configurations), which need to 
be factorized subsequently (i.e., approximated with some 
products of probabilities of single-site configurations). 
Eventually, each term contributing to the evolution of 
single-site probabilities will contain the probability that 
a chosen site is in an appropriate state (for a given pro- 
cess) and the rate of that process. In some processes, 
neighboring sites also give some contributions. Using the 
following notation for a probability of a site being in par- 
ticular state 

X\— empty 

x 2 — nutrient only 

X3 — tumor only 

x 4 — both nutrient and tumor 

one arrives at the following equations describing the time 
evolution of the probabilites (densities) Xi, i = 1, 2, 3, 4: 



~dt 

dx 2 

~~dt 

dx 3 

dt 

dx4 

~dt 
where / 



(1 - p)(x 3 - fxxXi) - pgxx 
pgxx - (1 -p)fx 2 x 4 

■ (1 - p)(x 4 - x 3 + fxix 4 ) - pg(l + A)x 3 
(l-p)(fx 2 x 4 -x 4 )+pg(l + A)x 3 , (9) 



1 — (x 3 +Xj) 2 



and = . , ; 



xi+x 2 " a l+A(a:3+a;4) ' 

Solving numerically the above set of equations, we de- 
termined the single-site probabilites and thus the tumor 
(xt — x 3 + x 4 ) and the nutrient (x n — x 2 + x 4 ) densities. 
To control the accuracy of the numerical procedure, one 
can check whether the obtained probabilites satisfy the 
obvious normalization condition x 4 + x 2 + x 3 + x 4 = 1. 



III. NUMERICAL RESULTS 
A. d=l 

First, we describe the results for the d — 1 version 
of the model. The already reported results of simula- 
tions show jl3[ that, for a sufficiently large p, the model 
remains in an active phase with Xt > 0, which termi- 
nates at a critical point p c (which depends on A). For 
p < p c , the steady state of the model is an absorbing state 
x t — and x n = 1 (tumor cells die out and the system 
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FIG. 1. (Color online) The time dependence of densities of 
tumor cells x t (t) calculated for d = 1, A = 3, and several 
values of p. At the critical point (p = 0.1764), the power-law 
decay Xt ~ t~ a is observed with the exponent a = 0.58(1). 
Forp = 0.1765, 0.1764, and 0.1763, the simulations were made 
for systems with N = 10 s sites, and for other values of p, 
systems of size N = 10 7 were used. 



gets filled with nutrients). At p = p c , the model under- 
goes a phase transition from an active into an absorbing 
phase and we expect that Xt is the corresponding order 
parameter. As it is already known, models with a single 
absorbing state are expected to belong to the so-called 
directed-percolation (DP) universality class In this 
universality class and for d = 1, the critical exponent 
«])P = 0.159(1) and the decay of the order parameter 
upon approaching the critical point is described by the 
exponent Pdp — 0.276(1). Numerical simulations show 
[13| that for A = (i.e., with no preference for placing 
the nutrients) the model belongs indeed to the DP uni- 
versality class. Much different behaviour, however, was 
reported for A > 0, where Monte Carlo simulations give 
(3 = 1.0(1) and a = 0.60(5). 

The numerical results reported in the present paper 
confirm the earlier findings concerning the critical be- 
haviour of this model. In particular, more extensive 
simulations (by an order of magnitude both in simula- 
tion time and system size N) yielded a more accurate 
estimation of the exponent a = 0.58(1) (FigfTJ). From 
the finite-size scaling of Xt(t) at criticality, we estimate 
z = 0.94(2) (Fig(2J). With a similar procedure, the expo- 
nent can be estimated from the off-critical steady- 
state data (Figj3]). The best collapse is obtained for 
z/|l = 1.75(5). 

While the exponent v\\ — 1.75(5) is quite close to 
its DP counterpart Mmp = 1.73(1), the exponent z = 
0.94(2) clearly differs from z DP = 1.58(1) Let us 
also notice that from the scaling relation (JT]) and our es- 
timation of a and U\\, we obtain /3 = 1.01(5) and such a 
value agrees with the previous estimation j3 — 1.0(1) [IH 
based on the steady-state values of the order parameter 
Xt (estimations of j3 based on the scaling relation (0 usu- 
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FIG. 2. (Color online) The rescaled time dependence of tu- 
mor cells density xt{t) calculated at the critical point p — 
p c — 0.1764 for several values of N. For z — 0.94, the curves 
for different values of TV collapse on a single curve. The in- 
set shows the same plot on a normal (non-logarithmic) scale. 
(a = 0.58) 
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FIG. 4. (Color online) The size dependence of steady-state 
tumor cells density xt (averaged only over surviving runs) 
calculated at the critical point of d = 1 and d — 2 versions of 
the model. From the slope of the data and relation (©, we 
estimate /3/v±. 
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FIG. 3. (Color online) The rescaled time dependence of tumor 
cells density Xt(t) calculated for several values of p. For vu = 
1.75(5), the curves for different values of p collapse on a single 
curve, (q = 0.58) 



ally are considered to be more accurate than those based 
on the steady-state data). Using the estimated value 
of /?, from the finite-size behaviour of the tumor den- 
sity at the critical point (©, we estimate v± = 1.80(5) 
(Fig©, which is very different from the Directed Perco- 
lation value v±dp = 1.097(1) 1]. Within the error bars, 
the obtained exponents z, v\\, and v± satisfy the relation 

The spreading-dynamics simulations confirm the 
power-law behaviour (© with 6 = 0.30(1) (Fig©, S = 
0.160(2) (Fig© and z' = 1.60(2) (Fig©. Such val- 
ues are very close to the Directed Percolation exponents 
e DP = 0.313(1), 8 DP = 0.159(1) and z' = 1.58(1) Q. 
Moreover, from the decay of P(t —> oo) ©, we estimate 
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FIG. 5. (Color online) The time dependence of the number 
of tumor cells Xt(t) calculated for d = 1, A = 3, and sev- 
eral values of p. The initial configuration contained a single 
tumor cell (Xt = 1) and all sites filled with nutrients. From 
numerical data for p = 0.1764, we estimate that <d — 0.30(1), 
which is very close to the d = 1 DP value Qdp = 0.313(1). 



P' = 0.277(2) (Fig©, which is in a good agreement with 
the DP value. One can also check that /?', 8, and vn 
satisfy the scaling relation (©. 

The obtained results are in our opinion surprising. 
While some exponents (mi, 9, 5, (3' , and z') are very close 
to their DP counterparts, some others (/3, a, v±, and z) 
are clearly different. Particularly puzzling is the fact that 
z 7^ z' . While from the finite-size scaling at criticality 
(Fig© we obtain z = 0.94(2), the spreading dynamics 
(Fig© leads basically to the DP value z' = 1.60(2). To 
our knowledge, in up to now examined models with a sin- 
gle absorbing state those exponents are equal (although 
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FIG. 6. (Color online) The time dependence of the tumor- 
survival probability P(t) calculated for d = 1, A = 3, and 
several values of p. The initial configuration contained a single 
tumor cell and all sites filled with nutrients. From numerical 
data for p = 0.1764, we estimate that S = 0.160(2), which is 
very close to the d = 1 DP value 8 DP = 0.159(1)0. 
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FIG. 8. (Color online) The critical decay of the ultimate 
survival probability P(t —¥ oo). Both in d = 1 and d — 2, 
the exponent j3' (estimated from the slope of the data) is 
consistent with the Directed Percolation values. 
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FIG. 7. (Color online) The time dependence of the mean 
square spreading R 2 (t) calculated for d = 1, A = 3, and 
several values of p. The initial configuration contained a single 
tumor cell and all sites filled with nutrients. From numerical 
data for p — 0.1764, we estimate that z' = 1.60(2), which is 
very close to the DP value z' DP — 1.58(1) [J. 



in some models, which, strictly speaking, lack an absorb- 
ing phase, they were also reported to be different [13] )• 

Using our estimations of critical exponents one can 
check that they satisfy the generalized hyperscaling rela- 
tion © but only when spreading exponent z' is replaced 
with z — 0.94(2) estimated from the finite-size scaling 
analysis. A further verification of scaling assumptions 
leading to relations (JS])-© would be needed to explain 
such a behaviour of our model. 



B. d=2 

Further numerical simulations showed that similar crit- 
ical behaviour is observed in the d = 2 version of our 
model. Indeed, the simulations showed that the decay of 
the steady-state tumor cell density Xt seems to be charac- 
terized by j3 = 1.0 (Fig|9]). Let us notice that a mean- field 
approximation (|9]) remains in a reasonably good agree- 
ment with Monte Carlo data despite a low dimensionality 
[d = 2) of the model. 

Moreover, from the scaling relation (3 = una and the 
known estimate for the d = 2 Directed Percolation [l| 
V\\ = 1.30(1), we obtain a = 0.77. The time-dependent 
simulations confirm that at criticality (p — p c — 0.08371) 
the decay of tumor cells follows the power law Xt(t) ~ t~ a 
with a as 0.77 (FigHUJ). 

From the finite-size scaling of the off-critical behaviour 
of Xt (t) , we conclude that the exponent v\ | is indeed very 
close to the expected DP value 1.3 (Fig|TT|). Moreover, 
from the finite-size behaviour of the tumor density at 
the critical point ([2]) and taking j3 — 1.0, we estimate 
v ± = 1.0(1) (Fig©. 

We also used the spreading method. The obtained re- 
sults, as in the d = 1 version, are consistent with d = 2 
Directed Percolation. In particular, O seems to be very 
close to 0.23 (FigOU and S as 0.45 (FigHSJ. The numeri- 
cal data for the mean square spreading are not presented 
here but they are also consistent with d = 2 Directed Per- 
colation. From the decay of P(t — > oo) (H|), we estimate 
P' = 0.56(2) (FigJS]), which again is in a good agreement 
with the DP value. All three exponents j3', 5, and v\\ are 
very close to the DP ones and they satisfy the scaling 
relation ([7]). 

Moreover, we estimated the exponent z from the col- 
lapse of the time dependent Xt(t) at criticality and for 
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FIG. 9. (Color online) The steady-state density of tumor cells FIG. 11. (Color online) The rescaled time dependence of tu- 

x t as a function of p calculated using Monte Carlo simulations mor cells density xt(t) calculated for d = 2, A = 3, and several 

and the mean-field approximation ((9J for d = 2, A = 3. The values of p. For v\\ = 1.30(5), the curves for different values 

inset shows that at the critical point p = p c = 0.08371 the of p collapse on a single curve, (a = 0.77) 
density x t decays linearly (/? = 1). 
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FIG. 10. (Color online) The time dependence of densities 
of tumor cells Xt(t) calculated for d = 2, A = 3, and several 
values of p. At the critical point (p = 0.08371), the power- law 
decay xt ~ t~ a is observed with the exponent a = 0.77(1). 
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FIG. 12. (Color online) The time dependence of the number 
of tumor cells Xt(t) calculated for d — 2, A = 3, and several 
values of p. The initial configuration contained a single tumor 
cell (X t — 1) and all sites filled with nutrients. The simula- 
tions show that at criticality (p — 0.08371) the exponent is 
very close to the DP value Qdp = 0.23(1). 



various values of N (Fig fH)) . Let us notice that our 
estimate z — 1-4(1) (finite-size scaling), together with 
a = 0.77(1), S = 0.45(1), and 6 = 0.23(1), satisfies 
the hyperscaling relation Such a value of z is also 
marginally consistent with the relation z = v\ i © . 

Finally, we would like to comment on the fact that 
both in the d = 1 and d = 2 versions of our model, the 
exponent /? most likely takes a unit value that is charac- 
teristic to infinite-dimensional (mean-field) models with 
absorbing state(s). In our opinion, it might be related 
to a roulette-wheel selection used to choose the site for 
placing a nutrient. This frequently used, e.g., in genetic 
algorithms, selection mechanism is in a sense nonlocal. 
For example, the removal of a tumor from a certain site 
changes the weight of this site, which affects the selection 



probabilites of all sites. This is because in the roulette- 
wheel method the selection probability of a certain site 
is given as a ratio of the weight of this site and the to- 
tal sum of all weights (of all sites). Thus a change at a 
certain site changes the total sum and that affects all the 
sites. Perhaps investigation of similar models but with 
different selection mechanisms would be needed to check 
whether the roulctte-wheel selection is indeed responsi- 
ble for such a critical behaviour. But even if the critical 
behaviour is changed due to the non-locality of roulette- 
wheel selection, some open questions still remain: (i) why 
only j3 acquires the mean-field value; (ii) why some oth- 
ers exponents (a, z, and v±) are changed but not to the 
mean-field values; (iii) and why some others (z/m, O, 8, 
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FIG. 13. (Color online) The time dependence of the tumor- 
survival probability P(t) calculated for d = 2, A = 3, and 
several values of p. The initial configuration contained a single 
tumor cell and all sites filled with nutrients. At criticality 
(p = 0.08371), the exponent S is very close to the DP value 
S DP = 0.45(1). [J]. 
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FIG. 14. (Color online) The rescaled time dependence of tu- 
mor cells density x t (t) calculated at the critical point p — 
p c = 0.08371 for d = 2, A = 3, and several values of N. 
For z = 1.4, the curves for different values of TV collapse on 
a single curve. The inset shows the same plot on a normal 
(non-logarithmic) scale, (a = 0.77) 



j3', and z 1 ) are intact and get the DP values. 



C. d=3 

As a last of our results, we present the steady state 
behaviour of our model in the d = 3 version. Qualita- 
tive dependence of x n and Xt on p (Fig fT5|) is similar to 
lower-dimensional versions. As expected, one can notice 
a much better agreement with mean-field approximation 
©• However, one can also observe that region of exis- 




FIG. 15. (Color online) The steady-state density of tumor 
cells xt and of nutrients x n as a function of p calculated using 
Monte Carlo simulations and the mean-field approximation 
© for d = 3, A = 3 and N = 100. 



tence of the absorbing phase (x n = 1, Xt = 0) shrinks. 
It means that even for low nutrient update rate p tumor 
cells survive in the system. Thus starving tumor to death 
is less effective in d = 3 systems but unfortunately this is 
the dimension that is most relevant from biological and 
clinical point of view. 



IV. CONCLUSIONS 

We presented a more detailed analysis of the critical 
behaviour of a recently introduced tumor growth model. 
Although the model has a single absorbing state, its criti- 
cal behaviour differs from the expected Directed Percola- 
tion universality class in an intriguing way In particular, 
only some of its critical exponents, namely j3, a, v>j_, and 
z, differ from their DP counterparts. Other exponents, 
namely v\\, ft 1 , and spreading exponents 0, (5, and z', 
most likely take the DP values. Moreover, the general- 
ized hyperscaling relation ([5]) seems to be satisified in our 
model, provided, however, that the spreading exponent 
z' is replaced with the dynamical exponent z. 

What is also interesting is the fact that in our model 
the exponent j3 seems to take the mean- field value /J = 1. 
In a speculative way, we attribute this fact to, in a sense, 
non-locality of the roulette- wheel selection, which is used 
in our model to choose a site to supply a nutrient. Even 
if it is so, it still remains unclear why such a non-locality 
affects only some of the exponents. 

Let us also mention that the fact that (3 — 1 exactly 
(if proven) might be used to precisely locate the critical 
point from the decay of the order parameter. With a 
precisely located critical point, one might expect that 
also other critical exponents, including DP ones, will be 
determined accurately. 

We hope that clarifying the behaviour of the present 
model might contribute to a better understanding of Di- 
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rected Percolation universality class, which is one of the 
most basic among nonequilibirum critical phenomena. 

Various models of tumor growth processes were used 
and shown to provide a valuable support to biological and 
clinical studies |22j . Often these models are very complex 



and their understanding is limited. Our approach uses 
a simplified model where more detailed understanding 
seems to be feasible. It would be interesting to examine 
more realistic extensions of our model that takes into 
account for example remodelling of vascular network that 
accompanies a tumor growth process (23j|. 
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